Zitterbewegung effect in spin-orbit coupled spin-1 ultracold atoms 
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The Zitterbewegung effect in spin-orbit coupled spin-1 cold atoms is investigated in the presence 
of the Zeeman field and a harmonic trap. It is shown that the Zeeman field and the harmonic 
trap have significant effect on the Zitterbewegung oscillatory behaviors. The external Zeeman field 
could suppress or enhance the Zitterbewegung amplitude and change the frequencies of oscillation. 
A much slowly damping Zitterbewegung oscillation can be achieved by adjusting both the linear 
and quadratic Zeeman field. Multi-frequency Zitterbewegung oscillation can be induced by the ap- 
plied Zeeman field. In the presence of the harmonic trap, the subpackets corresponding to different 
eigenenergies would always keep coherent, resulting in the persistent Zitterbewegung oscillations. 
The Zitterbewegung oscillation would display very complicated and irregular oscillation character- 
istics due to the coexistence of different frequencies of the Zitterbewegung oscillation. Numerical 
results show that, the Zitterbewegung effect is robust even in the presence of interaction between 
atoms. 

PACS numbers: 03.75.Mn, 67.85.Fg, 05.30.Jp, 



I. INTRODUCTION 



The Zitterbewegung (ZB) effect, which is character- 
ized by high frequency oscillations (trembling motion) for 
Dirac electrons was firstly predicted by Schrodinger [![. 
The Zitterbewegung effect is purely a relativistic phe- 
nomenon and originates from the interference between 
the positive and negative energy states of the electron. 
The experimental observation of the electron Zitterbe- 
wegung has not been realized due to its high frequency 
(the order of fruj = 2m e c 2 ~ lMeV) and small oscil- 
latory amplitude (the order of Compton wavelength of 
electron H/m e c ~ 10~ 12 m). However, it is shown that 
there exist Zitterbewegung- like effects in graphene 0, Q , 
superconductors photonic crystal @, single trapped 
ion and semiconductor quantum wells Hence, 
the Zitterbewegung-like effect has attracted great atten- 
tions recently, both theoretically and experimentally in 
various physics fields [HHL3j |. 



Most of these studies of the Zitterbewegung effect are 
based on the Dirac-like equation with spin-orbit coupling 
interaction. The spin-orbit coupling interaction plays an 
essential role in a lot of interesting phenomena, such as 
quantum spin Hall effects, topological insulator, exotic 
superconductivity or superfluidity, etc. In recent years, 
the spin-orbit coupling interaction in cold atoms have 
attracted great attentions both experimentally and the- 
oretically. The spin-orbit coupling in two component 
atoms have been created experimentally with Raman 
laser beams The proposal to realize spin-orbit 

coupling for three component atoms has been put for- 
ward in Ref. [18]. Using the so-called tetrapod setup 
scheme and two pair of counterpropagating laser beams, 
the spin-orbit coupling in spin-1 atoms could be real- 
ized in alkali-metal atoms. The resulting gauge potential 
is proportional to the projection of angular momentum 



operator of spin-1 atoms along xy plan. The ground 
states of the spin-orbit coupled BEC have been exten- 
sively studied theoretically. For example, in the case of 
strong coupling and the weak harmonic trap, the spin- 
orbit coupling would result in nontrivial ground state in 
BEC, such as plane wave phase or standing wave phase, 
which depends on interaction between atoms [19l - [21] . For 
the strong coupling and the strong harmonic trap, the 
half-quantum vortex phase or vortex lattice phase devel- 
ops [al HUM [Hi . 

The length and energy scales can be well controlled 
in cold atoms. Thus, it is possible to observe the 
Zitterbewegung-like oscillation in cold atom experiments. 
Recently, there were proposals to simulate the Zitterbe- 
wegung effect by using spin-orbit coupled ultracold atoms 
29 32], wherein most authors focus on the Zitterbewe- 
gung effect of two component atoms with spin-orbit cou- 
pling in free space. Thus, the quasi-momenta is a good 
quantum number and the Hamiltonian can be diagonal- 
ized within the momentum space. By using the Gaus- 
sian packet as the initial state, it is shown that the am- 
plitude of the Zitterbewegung oscillation deceases with 
time and the Zitterbewegung phenomenon has a tran- 
sient characteristic. It would be of interest to find ways 
to stabilize and manipulate the Zitterbewegung oscilla- 
tion. Zitterbewegung-like phenomena are fairly common 
characteristics for multi-level systems [33]. However, the 
multi-frequency Zitterbewegung oscillation in ultracold 
atomic system has not been investigated to the best of 
our knowledge. The effects which is induced by the Zee- 
man field and the trap on the Zitterbewegung oscillations 
need to be clarified. 

In the present paper, we investigate the characteris- 
tics of the multi-frequency oscillation induced by Zee- 
man field and the persistent oscillations in the presence 
of the external harmonic trap in the spin-1 ultra-cold 
atoms. In Sec. II, we introduce the general Hamiltonian 



2 



of spin-1 atoms with spin-orbit coupling, and then focus 
on the Zitterbewegung oscillatory characteristics in the 
Zeeman fields. In Sec. Ill, the Zitterbewegung oscilla- 
tory characteristics in harmonic trap are investigated. In 
Sec. IV, the effect of the interaction between atoms on 
the Zitterbewegung oscillation is considered. A summary 
is presented in Sec. V. 



spin-dependent interaction parameter, respectively, ao 
and 02 are the s-wave scattering lengths corresponding 
to the total spin of the two colliding bosons and 2, 
respectively. The spin-orbit coupling strength and the 
effective Zeeman parameters can be tuned by the wave 
length of Raman lasers and the detuning from Raman 
resonance, respectively [T3.fl5l|. 
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FIG. 1. (color online), (a) The amplitude of the Zitterbewe- 
gung oscillation as a function of the momentum fco under zero 
Zeeman field, (b) The Zitterbewegung oscillations under zero 
Zeeman field for various momentum fco- 



II. ZITTERBEWEGUNG EFFECT IN ZEEMAN 
FIELDS 

The two dimensional Hamiltonian of the spin-orbit 
coupled spin-1 atoms in the presence of the external Zee- 
man field and harmonic trap is 19] 



H = T+V t 



trap 



T -- 

Vtrap 



2m 



V so + V z 



Hi; 



muj 2 (x 2 + y 2 ) 



V SO = J dr^ + [ 1Px F x + iPyF y ]^ , 
V z = J df^+[ P F z + qF 2 }^, 
ffint= /df[|n 2 + |F 2 ], 



(1) 



where T, Vt rap , V S o, V z , H int are the kinetic energy, 
the harmonic potential, the spin-orbit coupling interac- 
tion, the effective Zeeman shift and the interaction be- 
tween atoms, respectively. * = (^i, ^o, ^-i) T , n = 
hi + n + h-i, F = V+FapVp, c = 47rft 2 (a + 2a 2 )/3m 
and C2 = 4:Trh 2 (a,2 — ao)/3 are the field operator, the 
number density operator, the angular momentum opera- 
tor, the spin-independent interaction parameter and the 



FIG. 2. (color online). Zitterbewegung oscillation in spin-1 
atoms under the linear Zeeman field. The Zitterbewegung 
oscillation under zero Zeeman field is plotted in green. The 
blue, red lines are Zitterbewegung oscillation under a linear 
Zeeman field p=4 and p=6, respectively. 
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The Pauli spin matrices for spin spin-1 atoms are given 
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Without the external harmonic trap and interaction, 
the Hamiltonian can be diagonalized in the momentum 
space as 



uix 
H' = [ u>2 
uj 3 



(2) 



where u>i = (p 2 +p 2 )/2m + ncosu — 6/3, 
Lo 2 = (p 2 x +p 2 y )/2m + ncos(u + 4tt/3) - 6/3, 
w 3 = (Px + Py)/^ m + ncos(u + 2n/3) — 6/3 are roots of 
eigen-equation Det(ujl — H) = 0. The corresponding 
eigenvectors are 



1 



~(pI +P 2 V ) + 2wi(> - q + Wi) 
\/2(p -q + Ui)(p x + ip y ) 

(Px + iPy) 2 



(3) 
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FIG. 3. (color online). Zitterbewegung oscillation under the 
quadratic Zeeman field. The blue and red lines are Zitter- 
bewegung oscillations corresponding to CJ23 and uj\2 for a 
quadratic Zeeman field parameter q=4. For the sake of com- 
parison, the green line corresponding to zero Zeeman field is 
also plotted. 



where rii is normalization coefficient, n = */— 4pi/3, u — 
arccos(- (Zl (-p 1 /3)- 3 / 2 /2)/3, Pl = c-b 2 /3, q x = d - 
bc/3 + 26 3 /27, b = -2q, c = -(pl + pg) + q 2 -p 2 ,d = 
(pl+p 2 )q, respectively. In general, they are the functions 
of the momentum (p x , p y ) and Zeeman parameters (p, q). 
The position operator in the Heisenberg picture is [33j 



f(t) = f(0) 



z, 



kd, 



(4) 



where V k = Qk = \a k )(a k \ are the group veloc- 

ity and the projection operator, respectively, k, I denote 
the energy branch u)x, u>2, ^3 and u>ki = ^k — &i is the 
eigen-energy difference. Z^a = iQk^§jf is the so-called 
Zitterbewegung amplitudes (see Ref. [331]). The first 
and second terms are constants in Eq. (4). The third 
term is the uniform motion and the fourth one corre- 
sponds to the Zitterbewegung oscillation. As shown in 
the fourth term, the position operator usually undergoes 
a multi-frequency oscillation. We will calculate the aver- 
age value of position by using an initial wave function of 
Gaussian density distribution: 



\g) = 



1 



x +y~ 



(5) 



The wave function expressed in momentum space is 

r [(p*-fc ) 2 +pJ] 1 n 1 (e) 



where S and kg are the width and the average momentum 
of the wave packet. The mean value of position is (f(t)} = 

(g\m\g)- " 




FIG. 4. (color online), (a) The energy spectrum under both 
the linear and quadratic Zeeman field (p=l, q=l). The green, 
red and black lines denote the three energy branches u>i, 012 
and W3, respectively. The black elliptic area is the region 
where the relative group velocity V12 approaches zero, (b) 
The Zitterbewegung oscillations for the momentum fco = 0.87. 
The red, black and green lines denote the Zitterbewegung 
oscillations with the frequency corresponding to uiy2 , CJ23 and 
cdi 3, respectively. 



From Eq. (4), we can find that the kinetic energy of 
Hamiltonian does not contribute to the Zitterbewegung 
oscillation. Therefore, during the discussion of the Zit- 
terbewegung in this section, the kinetic energy part is 
neglected. In this section, we take the reduced Plank 
constant H, the spin-orbit coupling strength 7 and the 
wave packet width S as independent fundamental units. 
The other derived physical quantities, such as time, mo- 
mentum, energy are measured by S/j, H/S and "fh/5, 
respectively. 

We explore the following four cases according to the 
applied Zeeman field. 

Case 1: The Zitterbewegung under zero Zeeman field 
(p = and q = 0). 

Considering the initial wave packet (5), only the ma- 
trix element of (y(t))i.i has contribution to the the oscil- 
lation part of y component of position operator. In the 
Heisenberg picture, it takes the form, 



y(th 



Px 



Pi + p; 



-cos(t<J : 



Px+P 



(7) 
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Substituting the initial state, we calculate the mean value 
of the oscillatory part of y 

(!/(*)) 



oo />oo 



oo -/ — oo 



dp— 

IT 



! Px cos(t^p 2 x + p 2) e -((Px-feo) 2 +pJ) 



/•OO 

= 2e- k »Re[ dpe- p2+up h{2kop)\, 
Jo 



(8) 



where p = yP 2 +P 2 and A (a;) are variable of integra- 
tion and the modified Bessel function of first order, re- 
spectively. In order to capture qualitative behaviors of 
Zitterbewegung oscillation, we assume » 1, and the 
asymptotic formula I\{x) m e x /y/2nx can be used f34| . 
The integral / °° dpe < --P 1+It P +2ko ^ /y/2ko~P can be evalu- 
ated using the method of steepest descents 



(»(*)) 



k 2 t 2 ; cos{k t- -), (9) 



where # = arcta^^j-) denotes a phase shift of Zitterbe- 
wegung oscillation. 

From Eq. (9), it is shown that in the case of fco ^> 1, 
the average position in the the direction perpendicular 
to the average momentum ko undergos a damping os- 
cillation with a single frequency. There are two fac- 
tors which result in the damping. One is the exponen- 
_* 2 

tially decreasing term e « which originates from the in- 
creasing spatial separation between the subpackets cor- 
responding to the higher and lower eigenenergy branches 
[35| . Because the overlap between the subpackets gets 
smaller and smaller with time, the amplitude of the os- 
cillation gets smaller and smaller. We can identify the 
non-vanishing relative group velocity between the sub- 
packets corresponding to different energy branches at the 
average momentum v — <9 Pa .wi2|fc = 1. It is anticipated 
that, when the relative group velocity gets smaller, the 
damping would be suppressed. The other one is the term 



V(fc 4 



) 4 , which results in a much slowly damping 



compared with the exponentially decreasing term. In the 
case of ko 3> 1, the exponentially decreasing term domi- 
nate the whole damping trend before the disappearing of 
the Zitterbewegung. Although there exists a phase shift 
8 during the damping process of Zitterbewegung oscilla- 
tion, we find that, in the case of ko ^> 1, the effect of 
the phase shift 8 is not obvious before the disappearing 
of Zitterbewegung. 

For the limit of fco <C 1, the modified Bessel function 
can be expanded as taylor series h(x) — § + fg + ^ + 
0(x 7 ). By keeping only the first linear term , we get 



(y(t)) « 2k e~ k ° / dpe~ p pcos(tp) 
Jo 

^koe-^il-tDi 1 -)), (10) 



where D(x) = ^ J °° e~* 2 l A s±n{xt)dt is the Dawson func- 
tion. 

We can see from Eq. (10) that there will be no inte- 
grate Zitterbewegung oscillation when the average mo- 
mentum is very small (see also panel (b) in Fig. 1). We 
notice that a similar result is found in 4 x 4 Luttiger 
Hamiltonian by Demikhovskii et al. 36] . 

From the above approximate formulas, the amplitude 
of Zitterbewegung oscillation y(t = 0) can be obtained for 
the two limit case. In the case of the average momentum 
ko <C 1, it is proportional to the average momentum. 
When ko 3> 1, the Zitterbewegung amplitude is inversely 
proportional to the average momentum. We numerically 
investigate the amplitude as a function of the average 
momentum in the panel (a) of Fig.l. In panel (a), there 
exists a maximum value of the amplitude on the curve. 
The maximum value occur at around fco ~ 1. 

Case 2: The Zitterbewegung under the linear Zeeman 
field (p ^ and q = 0). 

The oscillatory part of y(t) in the Hcisenberg picture 



is 



(y(t)h 



v 2 + P l +v 2 ^Pl+Vl+V 2 )- (11) 
rx * Vy ' tr 



In principle, one could be able to get a asymptotic re- 
sult though a similar calculations as the case of zero Zee- 
man field. However the resulting expression is so long 
and cumbersome that we could not get clear physical 
meaning from it. Inspired by the case of zero Zeeman 
field, we will give approximate expression to fit the data 
obtained from the exact numerical integral under condi- 
tions of fco 3> 1. The approximate formula is 



<„(t)> = (g\y(t)\g) 



2e~ k »Re[ I dp 



p e 



-h(2kop)] 



k 2 



K+p 2 (fc 4 + M^!)i 



(12) 



where v 12 = d Px L)i 2 \k 



_ kg 



is the relative group ve- 



locity at the average momentum and W12 = yfcjj + P 2 ^ s 
the energy difference between different energy branches 
at the average momentum. We can see from Eq. (12) 
that the amplitude is suppressed by the applied linear 
Zeeman field. With the increase of the linear Zeeman 
field, the decaying trend is suppressed. The reason is 
that with the increase of linear Zeeman field, the relative 
velocity between sub-packets corresponding to different 

energy branches at the average momentum V12 — , k ° 

V fc o+P 2 

gets smaller. Then, the sustained coherence between sub- 
packets lead to the suppression of decaying. For some 
specific parameters, we depict the Zitterbewegung oscil- 
lations with only the linear Zeeman field in Fig. 2. 

Case 3: The Zitterbewegung under the quadratic Zee- 
man field (p = and q ^ 0). 
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The oscillatory part y^i^of the position operator in 
the Heisenberg picture is 

,.s Px^23 COs(tWi 2 ) , P x 0Jl2 COs(tUJ23) M , 

where uj 23 = \{q + + P 2 y ) + q 2 ), w i2 = + 

^/4(p2+P?)+9 2 ), Wl3 = ^4(^+p2) +g 2. 

By using the similar approximation, the average posi- 
tion along y direction is approximated by 

(y(t)) = {g\y(t)\g) 

= Vl2 + 2/23 

v 2 t 2 

^2,3(^0) e 



2 ^13 (jfeg 



k 2 V 2 t 2 \ 
4 J 



-cos(wi 2 t) 



2u 



13 



(fco 4 



r cos(w 23 t), (14) 



where is the Zitterbewegung oscillation with fre- 
quency which is energy difference Wjj between the en- 
ergy branches i and j at the average momentum. V12 — 



2 k—kn 



V23 



9p^23 



are relative group velocities between different energy 
branches at the central momentum. 

In the presence of only the quadratic Zeeman filed, 
the Zitterbewegung oscillation split into two oscillations 
with different amplitudes and frequencies (see also Fig. 
3). The amplitude of the Zitterbewegung oscillation with 
higher frequency is smaller than the other one. Be- 
cause of the equal relative group velocities at the central 
momentum v\2 = V23, the two Zitterbewegung oscilla- 
tions have the same exponential decreasing factor and 
the decaying trends for the two Zitterbewegung oscilla- 
tions are the same on the whole. With the increase of the 
quadratic Zeeman field, the splitting of both the ampli- 
tudes and frequencies between the two Zitterbewegung 
oscillations gets more evident. 

Case 4: The Zitterbewegung under both the linear and 
quadratic Zeeman field (p 7^ and q 7^ 0). 

In the presence of both the linear and quadratic zee- 
man field, the analytical expression of the position op- 
erator in the Heisenberg picture is quite long and cum- 
bersome. We depict the Zitterbewegung oscillations with 
some specific parameters, e.g. p — 1, q = 1 and kg = 0.87 
in Fig. 4. We can see that there are usually three fre- 
quencies when both the linear and quadratic Zeeman 
field exist. We find that when the relative group veloc- 
ity i>ij = 9 Px u>ij I k between different eigen-states at the 
central momentum approaches zero, the mean position 
will damp much slowly. As shown in the bottom panel 
of Fig. 4, when the central momentum is chosen near 
the elliptic region and the relative velocity between sub- 
packets nearly vanish, the decaying trend of the Zitter- 
bewegung oscillation corresponding to wi 2 is suppressed 



fco 



greatly compared with the two others. The slowly de- 
caying amplitude may be favorable to the observation of 
Zitterbewegung in experiments. 



III. THE ZITTERBEWEGUNG EFFECT IN 
HARMONIC TRAP 

As shown in the above section, the amplitude of Zit- 
terbewegung oscillation in free space usually decays with 
time. In the trapped system, the oscillatory behaviors 
have unique characteristics compared with that in free 
space. Firstly, due to the confinement of trap, the atoms 
could not escape from the trap. Hence, the subpack- 
ets corresponding to different eigenenergies would always 
keep coherent, resulting in non-decaying Zitterbewegung 
oscillations. Secondly, due to the existence of infinite 
energy levels in the trap, the Zitterbewegung oscillation 
shows very complicated features. 

In this section, we adopt the natural units for har- 
monic oscillators. The length, mass and time are mea- 
sured by y/H/muj, m and 1/w respectively. The other 
physical quantities, such as, energy, momentum, velocity 
are measured by huj, y/hmuj and ^Jhui/m, respectively. 

Similar to that in Rcf. 



evolution operator as U 



33| . we decompose the time 



iHl 



,— iEtt 



where \l) is the eigenstates of Hamiltonian and Ei the 
corresponding cigenenergy. We generalize Eq. (3) to 
the trapped system and obtain the position operator in 
Heisenberg picture as 



r(t) = r*(0) 



k 



■ eiUMt Zk;l, 
k.l^k 



(15) 



where Z^-i = \k) (k\f\l) (l\ is the so-called Zitterbewegung 
amplitude operator, and uiki is energy difference between 
eigenenergies. To calculate the average value of posi- 
tion operator, the eigenstates and eigenenergies are re- 
quired. Before the numerical calculation of the eigen- 
equations, we will discuss the symmetries of Hamiltonian. 
The Hamiltonian have rotational symmetries along z ax- 
ial direction. Thus the total angular momenta along z 
direction J z = L z + F z is a good quantum number. We 
will label the eigenstates and energies with good quan- 
tum number j z = m. In the polar coordinates (p, 6), the 
eigen-function can be written in the following form 



-1)0 



2tt 

{m+i)fl 



(16) 



27V 



J 



with j z = 771. 

In addition to the rotational symmetry, there is time 
reversal symmetry in the absence of the Zeeman filed. 
The time reversal operator is expressed as T = UK with 
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FIG. 5. (color online). The Zitterbewegung oscillations with 
a Gaussian wave packet as an initial state in the trap for 
the wave packet width 5 = 1, the momentum ko = 2. The 
green, black and red lines correspond to the spin-orbit cou- 
pling strength 7 = 0.5, 2 and 8, respectively. 



U = exp(-iTrFy) and K the complex conjugate opera- 
tion. Its matrix form is 



T = 




(17) 



Then the time reversal state can be obtained as 
\tp- m (p,8)) — T\tjj m (p,9)} with degenerate eigenenergies 
E- m j = E m j. The eigenstates are doubly degenerate 
except for the states corresponding to m = 0. 

We will work with the basis of 2D harmonic oscillator, 
and the basis can also be introduced by two independent 
operators in terms of operators of a x 



37] 



1 



ad 



1 , 



y)i 

+ iCLy), 



(18) 



where a x r y \ is annihilation operators of the x(y) direction. 
The position and momentum operators can be expressed 
in terms of the operators in Eq.(18) and their correspond- 
ing adjoint operators. For example, the y component of 
the position is y — ^{ad — a J — a g + a+). The harmonic 



oscillator basis can be written as 

\Xnd,n g ~ 



(a+) n '(a+)*« 



|0,0), (19) 



where |0, 0) is the ground state of a 2D harmonic oscilla- 
tor. The 2D harmonic oscillator basis can be expressed 
in the coordinate space as 



{r\X 



rid,n g / 



Rn,m(p) 



2tt 



(20) 



where 



2(n!) 



"9 ■ 

\m\ 



m 

3 of' 



L 



Rn,m(p) 



is the associated Laguerre poly- 
is the binomial coefficient. 



(n+\m\) 

V" f~in—k (—x) 
n k=0^n+\m\ k\ 

nomial and C% = (n _ fc)!fc! 

From the expression of the position operators in terms 
of a g [d) an d their adjoint operators, there exist an impor- 
tant selection rule for Zitterbewegung amplitude opera- 
tors 

Z m ,i- m ',i' = \m,l)(m,l\r\m',l')(m',l'\ = 0, 

(m-ra'/il), (21) 

where I labels the eigenstates belonging to the same quan- 
tum m. It is shown that only for m = ml ± 1 , there 
exist no- vanishing Zitterbewegung amplitude. From the 
above equation, it is found that the second term in Eq. 



(15) vanishes, i.e. Z^ 



k;k 



0. 



One can express the Hamiltonian in the form of ma- 
trix by adopting the harmonic oscillator basis. The cor- 
responding eigenenergies and eigenstates can be obtained 
through direct numerical diagonalization (23j. 

Now we investigate the Zitterbewegung oscillatory 
characteristics in the trap. First of all, as shown by 
equation (15), when the the initial state is the superpo- 
sition of two eigenstates, there would be a an oscillation 
with one single frequency which is the difference of two 
eigen-energies. The amplitude of the oscillation does not 
decay. This is because that the components correspond- 
ing to two eigenenergies in the initial state always keep 
coherent in the trap as stated before. 

Next, we calculated the Zitterbewegung oscillations 
when an Gaussian wave packet with non- vanishing av- 
erage momentum is used as an initial state. The pro- 
jections of the initial state onto the harmonic basis are 
required in the calculation the average value of the po- 
sition. The detailed calculations of the projection are 
listed in the appendix. From Eq. (15), we know that the 
number of the frequencies of the oscillation involved in 
Zitterbewegung oscillation is usually arbitrary when the 
initial state is chose in a completely arbitrary fashion. 
Using the Gaussian wave packet as initial state, we find 
that there are indeed a lot of frequencies in the Zitterbe- 
wegung oscillation. 

For some specific parameters, the Zitterbewegung os- 
cillations in the trap are shown in Fig. 5. From Fig. 5, 
we can see that, for both the weak and strong spin-orbit 
coupling, the Zitterbewegung oscillations display the har- 
monic oscillator characteristics within short time interval 
(20 time unit in Fig. 5). This is because when the spin- 
orbit coupling approaches to zero, the energy-spectrum 
approach the harmonic oscillator spectrum. In the case 
of the strong spin-orbit limit, the eigenenergies can be ap- 
proximated by E m>n = [— 7 2 + 2n + l+m 2 /7 2 ]/2, which is 
analogous to the energy spectrum formula of two compo- 
nent atoms with strong spin-orbit coupling and harmonic 
trap [22l. [27l] . Compared with the two components atoms, 
the ground state is not degenerate. The energy spectrum 
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FIG. 6. (color online). The Zitterbewegung oscillations in the 
presence of both the harmonic trap and the Zeeman field with 
the wave packet width 5 = 1, the average momentum ko — 2 
and the spin-orbit coupling strength 7 = 0.5, respectively. 
The green, red and blue lines correspond to the Zitterbewe- 
gung oscillations for various Zeeman parameters. 



in strong spin-orbit coupling limit form landau level-like 
spectrum. Thus, in cases of both weak and strong spin- 
orbit coupling limit, the energy level spacing approaches 
hid. With the increase of the time, the out of phase oscil- 
lations begin to appear. Meanwhile, the oscillations with 
high frequency begin to manifest themselves (see the red 
line in Fig. 5). 

In the presence of the intermediate-strength spin-orbit 
coupling, the Zitterbewegung shows quite complicatedly 
and irregularly oscillatory characteristics. It is due to the 
complicated energy spectrum and coexistence of oscilla- 
tions with various frequencies in Zitterbewegung. 

We explore the Zitterbewegune oscillation in the pres- 
ence of both the Zeeman field and the harmonic trap. 
For the weak spin-orbit coupling, the top panel in Fig. 
6 shows that the Zitterbewegung effect is manifested 
as beat oscillation. The beat period gets larger and 
larger with the increase of the linear Zeeman filed. The 



FIG. 7. (color online). The Zitterbewegung oscillations in the 
presence of both the harmonic trap and the Zeeman field with 
the wave packet width 5=1, the average momentum ko — 2 
and the spin-orbit coupling strength 7 = 5, respectively. 



quadratic Zeeman field usually suppress the amplitude of 
the Zitterbewegung oscillation, as shown by the the mid- 
dle panel in Fig. 6. In the presence of both the linear and 
quadratic Zeeman fields, we can see that the beat period 
gets even much larger compared with the case where only 
the linear Zeeman field is included as shown the bottom 
panel in Fig. 6 (see the blue lines) . For strong spin-orbit 
coupling, the amplitude of the Zitterbewegune oscillation 
is usually smaller than that of weak spin-orbit coupling 
as shown in Fig. 7. At the same time the high oscillations 
is manifested more evidently. 



IV. EFFECT OF INTERACTION ON THE 
ZITTERBEWEGUNG OSCILLATION 

In this section, we numerically simulate the dynam- 
ical evolution of the spin-orbit coupled 23 Na atoms by 
solving the corresponding GP equation. Taking the ax- 
ial trap frequency uj z — 2tt x 800Hz, the transverse trap 
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FIG. 8. (color online). The Zitterbewegung oscillations in 
the presence of the both harmonic trap and the interaction 
between atoms with the wave packet width 5 = 1, the mo- 
mentum fco = 2, respectively. 

frequency w± = 2ir x 90Hz, the total number of atoms 
N= 10 4 and the experimental values of scattering length 
of 23 Na atoms ao = 46as, a 2 = 52as (0,3 is the Bohr 
radius) , the corresponding two dimensionless interaction 
parameters cq = 178.90, C2 = 7.16 (in natural units of 
harmonic oscillator). The effect of interaction on the 
Zitterbewegung oscillation in the trap is considered in 
Fig. 8. Comparing it with Fig. 5 (without interaction), 
the beat oscillation period get longer for weak spin-orbit 
coupling. For intermediate and strong spin-orbit cou- 
pling strength, the oscillatory amplitude for a long time is 
slightly suppressed. Even though the ZB effect originates 
from the single particle Hamiltonian, it is still robust in 
the presence of interaction between atoms. Therefore, 
the ZB oscillation is an universal phenomenon in the dy- 
namical evolution of the spin-orbit coupled atoms. 

The Zitterbewegune oscillation is accompanied by a 
damped oscillation of the probability density with time 
for each component (35j |. The oscillatory properties of 
position operator, such as the decaying trend, the os- 
cillatory frequency, e.g, can be reflected by the damped 
oscillation of the number of atoms. The oscillation of the 
number can be obtained through time-of-flight measure- 
ment [17] . The parameters involved in the Hamiltonian 
of atoms can be precisely controlled in cold atoms exper- 
iments. Taking 23 Na atoms for example, one can adopt 
the typical experiment parameters, such as the wave 
length of Raman laser A ~ 600nm, which corresponds 
to a spin-orbit strength 7 = y/2nh/ '(m\) rs 1.93cm/ s 
(m = 3.82 x 10~ 26 kg for 23 Na). One can prepare a gaus- 
sian packet of cold atoms and apply a series of time-of- 
flight measurement for various values of time. Then it is 



expected to obtain the decaying oscillation of the number 
of atoms. 



V. SUMMARY 

In summary, we have investigated the Zitterbewegung 
effect in spin-1 atoms in the presence of the Zeeman field 
and the external harmonic trap. It is shown that the Zit- 
terbewegung oscillations could be greatly affected by the 
external Zeeman field and the trap. The external Zee- 
man field could suppress or enhance the Zitterbewegung 
amplitude and change the frequencies of oscillation. In 
addition, multi-frequency oscillations could appear in the 
presence of the Zeeman field. Through adjusting both the 
linear and quadratic Zeeman field properly, we could ob- 
tain a much slowly damping Zitterbewegung oscillation. 
In the presence of harmonic trap, there would be the Zit- 
terbewegung oscillation without damping and with arbi- 
trary number of frequencies. The Zitterbewegung in the 
trap displays very complicated and irregular oscillation 
characteristics due to its complicated energy spectrum. 
It is shown that the Zitterbewegung oscillation is robust 
even in the presence of interaction between atoms. 

The present work should help us understand the roles 
of external Zeeman field and trap in the Zitterbewegung 
oscillations. The investigation of the effects induced by 
the external field and trap on Zitterbewegung oscillation 
opens up new possibilities for the manipulation and con- 
trol of the Zitterbewegung oscillation. We anticipate the 
Zitterbewegung effect be detected in the widely-studied 
cold atoms experiments with spin-orbit coupling interac- 
tions. 
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Appendix A: The projections of the Gaussian packet 
onto the harmonic oscillator basis 

In this appendix, we give the calculation of the pro- 
jections of the initial Gaussian wave packet onto the har- 
monic oscillator basis, which is required in the calculation 
of the average values of position operator. The initial 
wave packet reads: 



\g) = -JL=e-^ e ifc °*, (Al) 

where S denotes the width of the wave packet and fco is 
the average momentum of wave packet. 
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The projection onto the harmonic oscillator basis is 



bn,m\g) = 



2(_i)"jH / n \ 



(n + |m|)! 



+00 



p \rn\+l e -^ L \m\ {p 2 )JlmlMdp) 



(A2) 



where \<t> n . m ) = R n ^ m {p)S-j= i a 2 = |(1 + ^) and 

J| m |(a:) is the Bessel functions of order \m\. 

The associated Laguerre polynomial consists of various 
monomials. The integrals corresponding to the monomi- 
als can be obtained through an integral identity 



where T(z) is Euler's Gamma function and iFi(a; /3;x) 
is the generalized hypergeometric functions (38|. Com- 
bining Eq. (A2) with Eq. (A3), we can obtain the pro- 
jections by summing up all integral values corresponding 
to the monomials. The Zitterbewegung amplitude op- 
erators Zk-i in Eq. (15) are also calculated within the 
harmonic oscillator basis and the mean value of the po- 
sition operator can be obtained correspondingly. 



+ 00 



p»- 1 e -P 2 P 2 J v (ap)dp 



,V — jJL 



2 ' (_i)^ e 4F 1 F l (i—r. + i. v + i 



2p^T{u + I) 2p J 



— ) 
(A3) 
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